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ABSTRACT 

Space-based  object  detection  and  tracking  represents  a  fundamental  step  necessary  for  detailed  analysis  of  space 
objects.  Initial  observations  of  a  resident  space  object  (RSO)  may  result  from  careful  sensor  tasking  to  observe  an 
object  with  well  understood  dynamics,  or  measurements-of-opportunity  on  an  object  with  poorly  understood 
dynamics.  Dim  and  eccentric  objects  present  a  particular  challenge  which  requires  more  dynamic  use  of  imaging 
systems.  As  a  result  of  more  stressing  data  acquisition  strategies,  advanced  techniques  for  the  accurate  processing  of 
both  point  and  streaking  sources  are  needed.  This  paper  will  focus  on  two  key  methods  in  computer  vision  used  to 
determine  interest  points  within  imagery.  The  Harris  Comer  method  and  the  method  of  Phase  Congruency  can  be 
used  to  effectively  extract  static  and  streaking  point  sources,  and  to  indicate  when  apparent  motion  is  present  within 
an  observation.  The  geometric  inferences  which  can  be  made  from  the  resulting  detections  will  be  discussed, 
including  a  method  to  evaluate  the  localization  uncertainty  of  the  extracted  detections  which  is  based  on  the 
computation  of  the  Hessian  of  the  detector  response.  Finally,  a  technique  which  exploits  the  additional  information 
found  in  detected  streak  endpoints  to  provide  a  better  centroid  in  the  presence  of  curved  star  streaks  is  explained  and 
additional  applications  for  the  presented  techniques  are  discussed. 

1.  INTRODUCTION 

There  is  much  room  for  methods  of  computer  vision  and  image  processing  to  increase  the  utility  of  existing  space 
object  tracking  sensors.  Space  object  tracking  imagery  contains  rich  information,  not  all  of  which  is  being  exploited 
effectively.  As  a  result  of  image  processing  algorithms  being  written  to  handle  a  specific  type  of  acquisition 
strategy,  sensor  sites  may  be  limited  in  the  types  of  observations  and  resulting  data  products  which  may  be  reported 
from  the  data  they  collect.  There  have  been  several  papers  [1,2]  which  have  discussed  front  end  image  processing 
architectures  which  process  imagery  collected  in  sidereal  stare  [1]  and  rate  track  [2]  modes.  What  is  desired  is  a 
general  image  processing  architecture  capable  of  effectively  analyzing  imagery  collected  over  any  motion 
achievable  by  the  optical  system.  To  this  end  computer  vision  techniques  may  be  used  to  infer  general  rotational 
motion  of  an  optical  system  with  respect  to  a  star  field.  This  paper  seeks  to  motivate  the  use  of  additional  image 
features  to  perform  online  geometric  analysis  of  detected  signals,  enabling  when  possible  more  than  just  inferences 
of  detected  angles,  angle  rates,  and  brightness  of  objects  determined  by  centroids. 

The  geometric  analysis  introduced  in  this  paper  will  be  shown  to  be  effective  in  the  presence  of  stray  light,  or  non- 
uniform  scene  energy,  as  well  as  enable  on  the  fly  determination  of  sensor  rotational  motion  and  discrimination  of 
resident  space  objects.  The  proposed  processing  is  motivated  by  the  recent  proliferation  of  small  telescope 
technology  which  makes  new  sensor  architectures  possible  which  may  operate  in  a  dedicated  tracking  and 
characterization  mode  where  imagery  could  be  taken  at  video  rates.  The  rest  of  this  paper  is  organized  as  follows: 
section  2  will  describe  two  methods  of  signal  extracted  not  based  on  intensity  thresholds  but  on  characteristics  of 
locally  compact  signals.  Section  3  will  discuss  a  method  for  the  determination  of  centroid  and  endpoint  uncertainty 
using  the  Hessian  of  the  detector  output.  Section  4  presents  an  approach  to  determining  the  rotation  axis  and  rate  of 
a  slewing  imaging  system  which  may  be  used  to  simultaneously  identify  candidate  RSO’s  for  subsequent  tracking 
and  characterization.  Section  5  presents  a  method  for  mapping  image  data  to  a  new  representation  in  polar 
coordinates  centered  on  the  determined  axis  of  rotation.  In  this  representation,  a  match  filter  can  be  easily  derived 
for  the  identification  and  extraction  of  streaking  star  signals  which  may  be  both  dim  and  curved  in  the  original 
image. 
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2.  SIGNAL  EXTRACTION  BASED  ON  LOCAL  CURVATURE  AND  LOCAL  ENERGY 

There  is  a  large  body  of  work  dedicated  to  the  extraction  of  signals  from  imagery  for  the  purposes  of  space  object 
tracking.  A  key  assumption  in  many  methods  is  that  the  segmented  pixels  associated  with  a  given  object  may  be 
“centroided”  and  the  resulting  vector  will  be  the  direction  to  the  space  object.  For  unresolved  objects  the  center  of 
brightness,  or  intensity  weighed  centroid  is  used  to  derive  the  vector  direction  to  stars  and  RSOs.  As  objects  become 
partially  resolved  or  fully  resolved  however,  the  center  of  brightness  may  no  longer  represent  the  actual  direction  to 
the  center  of  mass  of  the  space  object.  Thus,  multiple  methods  of  the  extraction  of  information  of  detections  can  be 
used  to  infer  characteristics  of  observed  objects.  Two  methods  will  be  presented  here.  One  method  is  based  purely 
on  image  intensity  gradients  and  is  considered  a  standard  in  feature  point  extraction  for  object  tracking  in  computer 
vision.  The  second  method  is  based  on  frequency  domain  arguments  and  has  many  desirable  behaviors  for 
processing  space  imagery.  In  the  case  of  perfect  observation  of  a  static  source,  centroids  based  on  intensity,  local 
curvature,  and  local  energy  reside  in  the  same  place.  When  these  detection  locations  are  not  in  close  agreement,  it 
can  be  used  as  an  indicator  of  non-ideal  observation  conditions  which  may  cause  a  need  for  deeper  analysis  of  the 
detection. 

An  elegant  approach  to  comer  detection  was  presented  in  1988  by  Harris  and  Stephens  [3].  In  this  work  they  showed 
that  typical  local  minimum  of  the  auto-correlation  function  required  of  a  comer  or  a  “good  point  to  track”  is  easily 
determined  by  the  image  gradients.  Instead  of  comparing  windows  over  a  predefined  neighborhood  about  an  image 
pixel,  they  define  the  stmcture  tensor  over  a  point  neighborhood  and  introduced  a  Gaussian  weighting.  By 
evaluating  the  eigenvalues  of  the  stmcture  tensor  one  may  classify  pixels  as  belonging  to  a  local  edge  response, 
comer  response,  or  bland  region.  For  all  participating  pixel  members  of  a  Gaussian  weighted  neighborhood,  the 
weighted  sum  of  squared  differences  is  given  by: 

SSD(x,y)  =  EE  +  hV  +  i)  -  I (ij))2 

j  (i) 

Here  (i;  j)  are  the  vertical  and  horizontal  displacement  indices  of  the  pixel  neighbors.  The  weighting  function  w(i,  j) 
is  a  Gaussian  weighting  function  chosen  to  match  the  size  of  the  pixel  neighborhood.  The  displaced  pixel  I(x  +  j,  y  + 
i)  intensity  can  be  approximated  by  using  the  Taylor  Series  expansion  as: 

I{x+j,y  +  i)  -  I{ij)  +  Ix(ij)x  +  Iy{ij)y  (2) 

By  substitution  we  obtain: 

SSD(x,  y)  =  EE  w(i,j)(Ix(i,j)x  + Iy(ij)y)2 

■i  i  (3) 

The  expression  above  is  a  weighted  quadratic  function  of  the  image  partial  derivatives  and  (x,  y)  displacements  of 
the  pixel  neighbors.  This  can  be  represented  more  generally  in  matrix- vector  form  by: 

SSD(x,y)  =  [x  y]A[x  y]T  (4) 

The  eigenvalues  of  the  matrix  A  provide  the  proportions  of  the  principle  curvatures  of  the  local  autocorrelation 
function.  Simply  put  these  values  provide  a  direct  means  of  evaluating  the  local  curvature  of  image  data.  For  a 
comer  response,  both  eigenvalues  should  be  high;  this  implies  a  direct  relationship  to  the  minimum  eigenvalue 
detection  method  also  common  in  feature  extraction.  A  pixel  residing  on  an  edge  response  will  have  only  one  large 
eigenvalue.  Finally  a  pixel  located  in  a  bland  region  will  have  small  eigenvalues.  The  Harris  response  function  was 
developed  to  quickly  compare  the  eigenvalues  of  the  stmcture  tensor  and  determine  the  degree  to  which  the  data  was 
curved  and  therefore  represented  a  comer.  It  is  given  by: 
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i?  =  det(A)  -  kTr2(A)  =  aj3  -  k{a  +  p)  (5) 

Here,  k  is  an  empirically  determined  constant  used  to  scale  the  relative  size  of  the  eigenvalue  sum  compared  to  the 
product  of  the  eigenvalues.  By  applying  a  threshold  to  the  response  function  instead  of  directly  computing  and 
comparing  the  eigenvalues,  comer  feature  strength  can  be  determined  without  having  to  compute  the  more  costly 
square  roots  involved  in  the  eigenvalue  decomposition.  This  approach  works  for  images  of  any  resolution,  and  has 
seen  wide  application.  The  result  of  applying  the  Harris  comer  detector  to  an  image  of  a  star  field  is  shown  in 
figure  1 .  Notice  how  both  the  rate-tracked  RSOs  and  the  endpoints  of  the  streaking  stars  are  highlighted  in  the 
comemess  response.  The  addition  of  endpoint  information  to  standard  detection  or  centroiding  enables  the 
estimation  of  additional  parameters  of  the  observation.  In  this  case  the  right  ascension  and  declination  rates  of  the 
observation  can  be  measured  by  determining  the  streak  lengths  of  the  stars  as  well  as  the  ID’s  for  each  star.  An 
undesirable  trait  of  Harris  comer  localization  results  from  its  dependence  on  the  computation  of  image  gradients. 
Often  it  is  necessary  to  combine  a  smoothing  operation  with  the  computation  of  spatial  gradients  within  an  image. 
This  is  especially  required  in  noisy  imagery.  The  result  is  that  the  computed  peaks  can  be  slightly  biased  at  streak 
endpoints  where  the  signal  is  not  symmetric.  Care  must  also  be  taken  to  ensure  that  the  smoothing  kernel  scale  is 
carefully  selected.  Note  in  figure  1  that  the  RSO  signals  are  not  centrally  located  since  the  kernel  was  selected  to 
find  endpoints  of  the  streaking  stars.  For  imagery  where  significant  variation  in  star  and  RSO  signal  scales  can  be 
expected,  a  multi-scale  Harris  approach  is  required.  While  this  method  is  useful  in  determining  additional  features 
of  detected  objects,  it  is  not  immediately  amenable  to  direct  extraction  of  the  star  streaks  and  RSO  signals.  A  better 
method  for  direct  extraction  is  to  perform  segmentation  based  on  the  local  energy  model. 


Figure  1  -  Rate  track  image  of  4  GEO  objects  and  resulting  Harris  Comer  Response 


The  Phase  Congruency  model  for  feature  detection  [4]  is  not  based  on  image  gradients,  but  on  local  Fourier 
components  of  a  signal.  If  one  were  to  represent  a  location  in  an  image  by  a  local  Fourier  series  expansion,  that 
location  would  be  considered  perceptible  if  multiple  local  Fourier  components  were  in  phase.  Phase  Congruency  at 
any  angle  indicates  a  perceptible  feature,  and  the  angle  at  which  these  components  are  in  phase  can  be  used  to 
characterize  feature  type.  It  would  be  interesting  to  investigate  the  use  of  the  angle  of  congruency  as  a  method  for 
characterization  of  RSO’s  but  that  analysis  is  left  for  future  work.  This  characteristic  of  image  data  was  first 
analyzed  as  part  of  the  development  of  the  local  energy  model  which  has  been  used  by  Morrone  [5]  while  attempting 
to  explain  many  psychophysical  effects  of  human  perception.  The  practical  use  of  this  model  for  the  detection  of 
edges  and  comers  in  imagery  was  introduced  by  Kovesi  [4].  In  his  work  he  explains  that  one  may  interpret  the 
phase  congruency  in  any  location  of  the  image  as  the  ratio  of  the  radial  extent  of  the  path  taken  by  adding  the 
complex  valued  vectors  of  each  Fourier  component  in  the  phase  plane  head  to  tail,  to  the  overall  distance  traveled  by 
adding  the  vectors. 
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The  resulting  value  is  bounded  between  0  and  1  which  makes  interpreting  the  output  of  applying  Phase  Congruency 
to  and  image  convenient.  Additionally  it  is  straight  forward  to  apply  the  process  to  a  noisy  image  from  a  given 
sensor  and  determine  a  threshold  value  which  will  produce  a  desired  level  of  false  positives.  Figure  2  shows  the 
concept  of  local  frequency  components  being  in  phase  in  the  vicinity  of  a  feature  as  well  as  a  visualization  of  the 
geometric  interpretation  just  described. 


All  components  of  Square  Wave 
In  phase  at  Step  Edge 


Imaginary 


Figure  2  (left)  Multiple  Fourier  Components  of  a  periodic  signal  in  phase  at  the  step  edge,  (right)  Fourier 

components  plotted  tip  to  tail  in  the  polar  diagram 


In  very  much  the  same  way  as  comers  and  edges  were  determined  using  the  Harris  approach,  moments  may  be  taken 
of  the  Phase  Congruency  measure  and  used  to  indicate  the  same  stmcture.  One  may  build  the  matrix: 

r  pci  pc,  pcv 

■  [PCyPC,  PCl  J  g) 

In  this  case  the  minimum  and  maximum  singular  values  correspond  to  the  minimum  and  maximum  moments  of 
Phase  Congruency.  In  contrast  to  the  edge  and  comer  responses  of  the  Harris  method,  or  in  edge  and  comer 
detections  in  general,  the  comer  map  generated  by  this  algorithm  is  a  strict  subset  of  the  edge  map.  This  means  that 
extracted  edges  will  have  comers  associated  with  their  endpoints.  This  has  particular  use  in  this  application  because 
it  makes  the  output  immediately  capable  of  direct  extraction  of  static  or  trailing  point  sources.  The  correct 
extraction  of  stars  and  RSO’s  can  be  made  invariant  to  sensor  motion  which  is  very  desirable.  Figure  3  shows  the 
results  of  applying  the  Phase  Congruency  algorithm  to  the  same  image  as  was  used  in  the  Harris  detector  case.  Here 
a  threshold  is  set  on  the  minimum  moment  of  Phase  Congruency  in  order  to  distinguish  signal  from  noise.  Notice 
that  the  resulting  segmentation  mask  accurate  identifies  pixels  belonging  to  both  the  rate  tracked  RSOs  and  the 
streaking  stars.  Additionally  the  endpoints  of  the  streaking  stars  stand  out  as  additional  comers,  an  effect  shared  by 
the  Harris  detector.  In  contrast  however  one  may  use  this  method  to  perform  initial  signal  extraction,  and 
subsequent  geometric  analysis,  thus  providing  not  only  the  center  of  brightness  defined  by  the  intensity  of  the 
segmented  pixels  but  the  center  of  the  phase  congruency  weighting  and  any  other  additional  peaks  which  may  be 
important. 
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Figure  3  -  (Left)  Original  image.  (Center)  2nd  Moment  of  Phase  Congruency  indicating  comer  strength.  (Right) 
Resulting  segmentation  mask  derived  by  thresholding  the  2nd  Moment  value. 


Figure  4  demonstrates  the  robustness  of  this  method  to  stray  light.  Since  the  value  being  thresholded  is  a 
dimensionless  quantity  invariant  to  image  contrast,  the  output  performance  of  the  method  is  identical  even  though  a 
significant  ramp  has  been  added  to  the  background  noise  of  the  signal.  This  effect  is  very  desirable  for  developing 
an  image  processing  engine  that  is  capable  of  handling  varying  background  noise  and  stray  light.  In  the  next 
sections,  we  will  demonstrate  how  the  output  of  these  methods  can  be  used  to  define  detection  uncertainty  and  to 
infer  both  sensor  and  object  motion  from  imagery  alone. 


Figure  4  -  (Left)  Original  image  with  addtion  of  discontinuous  non-uniform  noise  floor.  (Right)  Result  of  Phase 
Congruency  computation.  The  structure  of  the  output  is  invariant  to  the  presence  of  the  background  variation. 


3.  CENTROID  AND  ENDPOINT  UNCERTAINTY  DETERMINATION 

The  direct  computation  of  uncertainty  in  point  based  features  is  well  known  in  computer  vision.  While  many 
methods  which  use  principle  component  analysis  (PCA)  or  information  contained  in  the  eigenvalues  computed  from 
the  positions  and  values  of  segmented  pixels  exist,  a  direct  relationship  between  the  local  Hessian  of  the  detector 
response  function  and  the  localization  uncertainty  has  been  derived  and  effectively  used  for  many  feature  extraction 
and  tracking  applications  [6].  The  hessian  is  the  matrix  of  second  order  partial  derivatives  which  may  be  computed 
for  every  pixel  of  any  image.  For  example  localization  uncertainty  is  expressed  as  the  inverse  of  the  local  hessian 
matrix  of  the  Phase  Congruency  response: 
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Here,  1  and  1’  may  take  both  values  x  and  y  in  order  to  fill  out  the  four  partial  derivatives  which  comprise  the  matrix. 
Figure  5  demonstrates  the  computation  of  localization  uncertainty  at  every  pixel  of  the  Phase  Congruency  algorithm 
response  computed  in  the  vicinity  of  a  streak  (highlighted  by  the  box  in  figure  4).  The  minimum  localization 
uncertainty  coincides  with  the  peak  values  detected  at  the  streak  endpoints. 


Figure  5  -  (Left)  Local  neighborhood  of  Original  image.  (Center)  Local  Phase  Congruency  Computation.  (Right) 
Uncertainty  calculated  for  each  pixel  of  sufficient  response.  Minimum  covariance  is  located  at  the  detection  peaks. 


In  stark  contrast,  one  may  use  a  match  filter  which  assumes  a-priori  knowledge  and  evaluate  the  a-posteriori 
localization  uncertainty  by  computing  the  Hessian  of  the  match  filter  response  as  shown  in  figure  6. 


Figure  6  -  (Left)  Local  neighborhood  of  Original  image.  (Center)  Derived  Match  Filter  (Right)  Uncertainty 
calculated  for  each  pixel  of  sufficient  response.  Minimum  covariance  is  located  at  the  detection  peaks. 


Without  prior  information  detection  uncertainty  and  knowledge  of  the  signal  geometry  in  an  image  is  coupled.  An 
iterative  procedure  can  be  easily  designed  which  progresses  from  zero  a-priori  detection,  to  maximum  a-posteriori 
signal  extraction  by  minimizing  the  localization  uncertainty  and  maximizing  the  knowledge  of  the  observed 
trajectory  within  the  frame.  The  next  sections  will  show  how  this  information  can  be  used  to  effectively 
discriminate  RSOs  from  stars,  and  to  determine  important  parameters  describing  an  image  taken  from  a  slewing 
platform. 


4.  ONLINE  DETERMINATION  OF  FIELD  ROTATION  AND  DISCRIMINATION  OF  RSO’s 

With  a  set  of  endpoints  for  all  signals  detected  within  an  image,  we  are  able  to  infer  additional  information  about  the 
motion  of  the  sensor  as  well  as  the  particular  behavior  of  individual  signals.  To  achieve  this,  we  must  first  group  the 
endpoints  according  to  occurrence  in  time.  Note  that  there  is  an  ambiguity  in  which  endpoints  occurred  at  the  initial 
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time  and  which  correspond  to  the  final  time,  however  we  are  able  to  determine  which  endpoints  likely  occurred 
simultaneously.  This  solution  takes  a  geometric  approach  in  the  image  plane  to  determine  the  probable  intersection 
of  the  axis  of  rotation  with  the  image  plane,  hereafter  referred  to  as  the  center  of  rotation. 

Assuming  that  each  signal  corresponds  only  two  endpoints,  we  begin  by  computing  the  midpoint  of  the  endpoint 
pair,  ( xi)C,yi)C ),  where  i  =  1,2, ... ,  N  designates  the  signal  index.  With  this  point  and  the  slope  of  the  line  connecting 
the  endpoint  pair,  mh  we  are  able  to  write  the  equation  of  perpendicular  line  passing  through  the  midpoint  as 

y  -  Vi,c  =  —  (x  -  xi  c) 

m.  (9) 

Rearranging, 


1 

S'Lc  “r  VLc 

rrti 


(10) 


The  intersection  of  the  perpendicular  lines  for  each  streak  corresponds  to  the  location  of  the  center  of  rotation.  We 
are  able  to  build  a  least-squares  problem  of  the  form  AX  =  b, 
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(11) 


Where  X  E  R2  denotes  the  center  of  rotation.  The  least-squares  estimate  of  the  center  of  rotation,  then,  is  given  by 

X  =  (ATA)~1ATb  (12) 

To  correctly  group  simultaneously  occurring  endpoints,  we  compute  the  vector  from  each  midpoint  to  each  of  its 
corresponding  endpoints  and  consider  the  vector  from  the  center  of  rotation  to  each  midpoint,  as  seen  in  Figure  7. 


Figure  7  -  Vector  cross  product  operation  for  endpoint  grouping 


The  correct  endpoint  grouping  follows  from  a  cross  product  operation  between  Vc  and  the  endpoint  vectors,  V1  and 
V2 ,  for  each  streak.  A  sign  difference  in  third  component  of  the  result  is  an  indicator  of  endpoint  grouping,  i.e. 


[V«  x  V,]3  =  -\vc  X  Vj]3 


(13) 


Note  that  in  order  for  the  cross-product  to  be  defined,  all  three  vectors  must  be  augmented  with  an  arbitrary  third 
component.  For  scenarios  where  the  streaks  are  nearly  straight,  the  term/lM  approaches  singularity.  In  this  case, 
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however,  the  grouping  of  endpoints  is  fairly  trivial.  Note  that  any  center  of  rotation  estimate  beyond  the  bounds  of 
the  image  and  along  the  perpendicular  line  of  any  streak  will  provide  the  correct  endpoint  groupings  with  the  vector 
cross-product  procedure. 


The  correct  endpoint  grouping  now  allows  an  estimate  of  the  field  rotation  from  a  single  frame.  By  computing  the 
body-fixed  unit  vector  corresponding  to  each  endpoint,  we  are  able  to  formulate  Wahba’s  problem, 


(14) 


where  bt  and  rt  are  two  sets  of  unit  vectors,  A  describes  the  rotation  between  them,  and  at  is  a  set  of  weightings  [7]. 
Many  solutions  of  Wahba’s  problem  exist  which  will  provide  an  estimate  of  the  rotation  between  the  two  endpoint 
groups  [8,9].  This  solution  additionally  provides  a  method  of  detecting  RSOs.  Computing  the  residuals  of  the 
rotation  estimate  as 


e-i  = 


(15) 


yields  an  indicator  of  non-star  behavior.  Note  that  this  works  under  the  assumption  that  the  field  contains  mostly  star 
streaks  and  only  a  small  number  of  RSOs.  In  scenarios  with  closely-spaced  objects  and  a  small  field  of  view,  the 
rotational  estimate  may  become  biased  by  the  RSO  detections  such  that  it  becomes  unreliable. 


5.  STREAK  PROCESSING  IN  POLAR  COORDINATES 

Knowledge  of  the  sensor  rotation  allows  us  to  perform  a  more  thorough  analysis  of  the  imagery  to  discover 
additional  signals.  In  general,  it  is  possible  to  extract  dim  signals  from  an  image  with  knowledge  of  the  appearance 
of  the  signal  uncorrupted  by  noise.  This  template  matching  procedure  is  referred  to  as  a  matched  filter  [10].  From 
the  reference  frame  of  a  rotating  sensor,  however,  characterization  of  the  expected  size  and  shape  of  a  star  signal 
becomes  more  complicated.  In  the  case  of  Figure  8,  we  see  that  the  size  and  shape  of  a  star  streak  is  dependent  on  its 
location  within  the  image.  One  could  perform  and  adaptive  matching  process,  wherein  the  algorithm  generates  the 
expected  signal  for  each  location  in  the  image,  however  this  would  be  a  computationally  intensive  algorithm. 


Figure  8  -  Streaked  star  image  with  an  RSO 


Alternatively,  we  are  able  to  utilize  the  center  of  rotation  estimate,  derived  in  the  previous  section,  to  simplify  the 
pattern  matching  procedure.  By  re-sampling  the  image  in  polar  coordinates  about  the  center  of  rotation,  we  are  able 
to  produce  an  image  where  star  signals  are  equal  length,  and  that  length  corresponds  directly  to  the  angle  of  rotation 


DISTRIBUTION  STATEMENT  A:  Unlimited  Distribution 


DISTRIBUTION  STATEMENT  A:  Unlimited  Distribution 


of  the  sensor.  We  define  a  polar  coordinate  system  with  the  center  of  rotation,  X,  as  the  origin.  The  relation  between 
the  (x,y)  image  coordinates  and  (r,  6 )  coordinates  is 

x  —  rcos((9)  ■+-  X\ 

y  —  rsind  Xi  (16) 

Thus  a  grid  of  (r,  6 )  coordinates  can  be  converted  to  (x,y)  coordinates,  enabling  the  user  to  re-sample  the  original 
image  in  the  polar  coordinate  frame.  The  image  from  figure  8  is  seen  again  in  figure  9,  re-sampled  in  the  polar 
coordinate  frame  about  the  center  of  rotation.  Note  that  interpolated  values  of  regions  beyond  the  image  boundaries 
are  undefined.  These  should  be  indicated  as  such  and  removed  from  consideration  in  future  processing.  Here,  out-of- 
bounds  regions  are  marked  as  NaN  (not  a  number),  and  can  be  seen  in  the  lower  portion  of  the  image  as  dark  areas. 


Choice  of  gridding  parameters  depends  on  the  location  of  the  center  of  rotation  with  respect  to  the  image  plane.  The 
range  of  r  and  6  should  be  chosen  such  that  the  original  image  is  completely  sampled  with  as  little  out-of-bounds 
interpolation  as  possible.  It  is  also  worth  noting  that,  in  regions  near  the  center  of  rotation,  this  re-sampling  method 
produces  artifacts  from  the  image  noise.  This  effect  can  be  seen  in  the  upper  portion  of  figure  9.  For  small  radii 
about  the  center  of  rotation  the  image  is  sampled  too  densely,  leading  to  apparent  correlations  in  the  re-sampled 
image.  This  is  typically  only  a  problem  when  the  center  of  rotation  is  within  the  boundaries  of  the  image,  however  it 
can  also  arise  in  any  case  if  the  sampling  density  of  (r,  0)  is  high. 

The  re-sampled  image  allows  for  additional  analysis  which  is  not  possible  in  the  original.  For  example,  the  apparent 
motion  of  the  RSO  is  immediately  apparent  relative  to  the  stellar  background.  Since  all  of  the  star  streaks  in  the 
polar  frame  are  straight  with  identical  lengths,  any  streak  which  which  deviates  from  this  indicates  some  apparent 
motion.  We  can  take  an  additional  step  in  this  direction  by  applying  a  matched  filter  to  the  re-sampled  image.  By 
estimating  the  streak  length  from  the  re-sampled  image  (or  inferring  it  from  an  on-board  gyro),  we  are  able  to  create 
a  template  for  any  star  in  the  field  of  view.  Correlating  this  template  with  the  image  produces  the  result  seen  in  right 
of  figure  9.  Note  that  we  find  a  strong  response  for  all  star  signals  without  the  need  for  an  adaptive  matched  filter. 


6.  GEODETICA:  GENERALIZED  SSA  IMAGE  ANALYSIS  ENGINE 

The  Generalized  Electro-Optical  Detect  Track  ID  and  Characterize  Application  (GEODETICA)  is  a  data  reduction 
pipeline  which  combines  advanced  computer  vision  and  image  processing  algorithms  with  online  orbit 
determination  and  photometric  analysis.  The  Matlab  based  software  is  available  for  testing  and  has  been  used  to 
process  simulated  and  real  space  object  tracking  datasets.  In  addition  to  the  techniques  presented  in  this  paper, 
frame  to  frame  association  is  performed  using  a  Bayesian  association  method  and  Kalman  filter  for  tracking  of  RSO 
and  star  positions  as  well  as  an  online  orbit  determination  approach  which  uses  a  particle  filter  based  on  the 
constrained  admissible  region.  A  key  characteristic  of  the  image  processing  pipeline  is  the  use  of  multiple  feedback 
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mechanisms  between  the  image  processing  and  the  estimators  supporting  the  orbit  determination  and 
characterization  functions.  Additionally,  this  software  has  been  wrapped  with  the  Astrometry.net  star  ID  process 
which  enables  effective  determination  of  key  optical  parameters,  distortion,  and  inertial  pointing  for  a  wide  array  of 
ground  and  space  based  sensors. 


7.  CONCLUSIONS  AND  FUTURE  WORK 

This  paper  presented  two  key  methods  from  computer  vision  used  to  determine  interest  points  within  imagery.  The 
Harris  Comer  method  was  shown  to  provide  a  means  to  detect  streak  endpoints,  but  suffers  from  bias  problems  due 
to  the  smoothing  necessary  to  localize  the  peaks.  The  method  of  Phase  Congruency  was  shown  to  be  quite 
promising  in  providing  a  means  for  direct  space  object  segmentation  and  geometric  analysis.  Geometric  inferences 
which  can  be  made  from  extracted  signals  were  discussed,  including  a  method  to  evaluate  the  localization 
uncertainty  of  the  extracted  detections  which  is  based  on  the  computation  of  the  Hessian  of  the  detector  response. 
Finally,  a  technique  which  exploits  the  additional  information  found  in  detected  streak  endpoints  to  provide  a  better 
centroid  in  the  presence  of  curved  star  streaks  is  explained  and  additional  applications  for  the  presented  techniques 
are  discussed.  By  enabling  the  match  filter  in  polar  coordinates  it  is  possible  to  distinguish  very  low  SNR  stars  from 
RSO’s  in  a  single  frame  given  that  streaking  motion  is  apparent  over  the  exposure  time.  Future  work  will  focus  on 
the  development  and  maturation  of  GEODETIC  A,  an  image  processing  pipeline  developed  at  AFRL/RV  which 
combines  in  a  unique  way,  front  end  image  processing  with  orbit  determination  and  photometric  analysis.  It  is  the 
vision  of  the  authors  that  this  Matlab  toolbox  may  be  used  to  enable  robust  interpretation  of  imagery  taken  from 
ground  and  space  sensors  that  are  in  various  operating  modes,  and  that  the  list  of  available  options  for  sensor  CON- 
OPS  grows  significantly  by  having  a  generalize  framework  which  combines  image  processing  with  celestial 
mechanics  and  sensor  attitude  mechanics.  The  techniques  presented  here  a  just  a  few  of  the  low-level  processing 
algorithms  which  are  used  as  building  blocks  to  enable  efficient  interpretation  of  space  object  tracking  data. 
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